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Transition state trajectory stability determines barrier crossing rates in 
chemical reactions induced by time-dependent oscillating fields 
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When a chemical reaction is driven by an external field, the transition state that the system must pass 
through as it changes from reactant to product—for example, an energy barrier—becomes time-dependent. 
We show that for periodic forcing the rate of barrier crossing can be determined through stability analysis 
of the non-autonomous transition state. Specifically, strong agreement is observed between the difference in 
the Floquet exponents describing stability of the transition state trajectory, which defines a recrossing-free 
dividing surface [G. T. Craven, T. Bartsch, and R. Hernandez, Phys. Rev. E 89, 040801 (R) (2014)], and the 
rates calculated by simulation of ensembles of trajectories. This result opens the possibility to extract rates 
directly from the intrinsic stability of the transition state, even when it is time-dependent, without requiring 
a numerically-expensive simulation of the long-time dynamics of a large ensemble of trajectories. 


Controlling the rate at which reactants transform to 
products, either to accelerate a chemical process or to 
bias a reaction toward a certain pathway, is fundamen¬ 
tal to chemical physics. Such kinetic control can be 
achieved through forcing from an external field, leading 
to emergent behavior in molecular structure assembly^^— 
organic synthesis^, ultracold chemical reactions^ and sin¬ 
gle molecule spectroscopyii In these processes, reaction 
rates are typically obtained through transition state the¬ 
ory (TST)4“— There are two major obstacles to the im¬ 
plementation of TST. First, reactive trajectories must be 
identified and, second, the flux of these reactive trajec¬ 
tories though a phase space dividing surface (DS) must 
be calculated. If this DS is recrossed by reactive trajec¬ 
tories, TST overestimates the rate. Only in cases where 
this DS is recrossing-free is TST formally exact. 

In autonomous systems, the optimal DS is deter¬ 
mined by a normally hyperbolic invariant manifold 
(NHIM)i ^°d^~ — The study of NHIMs is the principle fo¬ 
cus of modern reaction dynamics in so far as knowledge of 
their geometry inherently contains the determining char¬ 
acteristics of the reaction. However, even when a re¬ 
crossing free DS can be found, a rate calculation can be 
intractable, especially for systems with many degrees of 
freedom, as large numbers of trajectories must be inte¬ 
grated to yield statistically relevant results. 

When a reaction is subjected to a time varying exter¬ 
nal force, the geometric structures of TST are known 
to exist in several cases, though they become time 
dependent For chemical reactions that are induced 
solely by by an external field, the coupling of the field 
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with the reacting molecule’s dipole moment can acceler¬ 
ate the reaction rate^^ even for systems that dissipate 
energy through a spontaneous emission processi^ 

An example of a molecular process where an external 
force influences the transition state geometry, and thus 
reaction rates, is the photoinduced isomerization between 
cis and trans stilbene (Ph-C=C-Ph)42r— Its unimolecu- 
lar reaction path can be parameterized through the tor¬ 
sion angle of the C=C double bond. Changing the ener¬ 
getics along this path through photoinduction alters the 
isomerization reaction rate. 

We show here that when a chemical reaction is peri¬ 
odically forced by an external field (such as a laser), the 
reaction rates are determined directly by the stability of 
the transition state. We calculate the reaction rate of a 
model system by simulating large ensembles of trajecto¬ 
ries and compare this result with the rate predicted by 
Floquet analysis of the transition state trajectory. Cor¬ 
responding to the “chemical method” where the reac¬ 
tant concentration is followed as a function of time^ 
we obtain reaction rates from the decay of a given ini¬ 
tial distribution. These rates are well-defined because 
the decay is exponential when averaged over a period of 
the driving and independent of the choice of distribu¬ 
tion. A major result of this work is that the rates can be 
obtained from a Floquet analysis of the transition state 
trajectory, an unstable periodic orbit close to the bar¬ 
rier top. This agreement suggests that chemical reaction 
rates can be extracted directly from the transition state 
without knowledge of the dynamics of the reactive pop¬ 
ulation. This general result could have been anticipated 
from the known connection between the stability of pe¬ 
riodic orbits of Hamiltonian systems and rates but 
is here established even in the case of driven systems. 

To model barrier crossings in chemical reactions driven 
by a time-dependent external field E{t) we consider a 
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FIG. 1. The time evolution of x{t) for a swarm of trajectories 
following the Eq. © with E(t) = asm{nt + tp) are shown in 
black (below). The potential surface of Eq. (IT|) is shown above 
with a contour plot shown below. The BT and TS trajectories 
are shown in dashed and solid white, respectively. Time is 
shown in units of r = Q,tj2'K + 3/4. The initial velocities are 
sampled from qb- Parameters are e = l,fl = 3,7 = 4, and 

ip = 0. 


particle of unit mass with an initial position xq on the 
reactant side of a moving energy barrier. The chosen 
barrier is a quartic potential of the form 

U{x) = -^uj^(x - E(t)f - }e(x - E(t))^, (1) 

which leads to the equations of motion 

X = V, 

( 2 ) 

V = -'yv + uj1{x - E(t)) + t{x - E{t))^, 

where 7 is a dissipative emission parameter, Wb is the bar¬ 
rier frequency, and e is an anharmonic coefficient. The 
anharmonic coefficient is restricted to values e > 0 such 
that there is a single maximum in the potential located 
at the barrier top (BT). The time dependent, instanta¬ 
neous position of the BT is specified by E{t). Figure [T] 
shows the time evolution of x{t) for an ensemble of tra¬ 
jectories following Eq. ©• Each trajectory either crosses 
the energy barrier forming product or remains on the re¬ 
actant side, never surmounting the barrier. The normal¬ 
ized flux of reactive trajectories through the phase-space 
bottleneck —the TS— is the reaction rate4 

Every realization of the forcing E{t) has a special tra¬ 
jectory imbedded in the dynamics I© that remains close 
to the BT for all time. This bounded trajectory, termed 
the transition state (TS) trajectory^^— will never de¬ 
scend into the product or reactant regions.— As illus¬ 
trated in Fig. [U the TS trajectory does not follow the 
time evolution of the energetic maximum given by the 
BT. It is instead a specific trajectory that responds to 
motion of the BT in such a way that it remains bounded 
for all time. When E(t) is a periodic function with period 
T such that E{t) = E{t + T) for all t, the resulting TS 



FIG. 2. Phase space plots for a swarm of trajectories following 
Eq. ([2} with E(t) = o sin(nt-|- 95 ). The initial position for ev¬ 
ery trajectory, xq, is shown as a dashed black line. Reactive 
trajectories are colored in blue and nonreactive trajectories 
are colored in orange with respective basins separated by 
(solid red). The TS trajectory is shown in black. The 
critical velocity is indicated by a red circle at the intersec¬ 
tion of the dashed red line and xq. The initial velocities are 
sampled from qb- Parameters are = 5, 7 = 2, and ip = 0. 


trajectory is a periodic orbit (PO) with the same period 

T. 

Attached to the TS trajectory are stable and unstable 
manifolds. The stable manifold intersects a line of initial 
conditions x = xq aX a critical velocity I/I— A tra¬ 
jectory will surmount the energy barrier, moving from 
the reactant state to the product state, if vq > I/I. If 
vq < I/I, the trajectory is nonreactive. The extension 
of this point to all values of Xg creates a critical curve 
V/I, which is a time-invariant phase space separatrix as 
illustrated in Fig. [21 Knowledge of I/I allows the identifi¬ 
cation of reactive trajectories from initial conditions, but 
it does not contain direct dynamical information such as 
the reaction rates themselves. 

To calculate rates, the TST methodology is concerned 
with creating a DS that is crossed once and only once by 
reactive trajectories and then evaluating the flux through 
that surface. For the case when Vj is known exactly, 
the no-recrossings criterion is satisfied and TST gives the 
formally exact reaction rate. In practice, large numbers 
of trajectories are generated and the flux is calculated 
through brute force. To construct a recrossing-free DS 
we will use a time-dependent DS that is located at the 
instantaneous position of the TS trajectory. As shown 
previously by uSf^^ the conhguration space projection of 
the TS trajectory is free of recrossings. 

For the case of a harmonic barrier (e = 0), Eq. © 
can be solved analytically with eigenvalues Au,s = 

— ^ ^7 ± y '72 -p corresponding to the unstable and 
stable manifolds, respectively. The TS trajectory is given 
in Refs. EH and El as 
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in terms of the S functionals^!^ 


ST[fJ-,g]t] = 


/ oo 

g{T) — t)) dr : Re^ > 0, 

+ / ^(t) exp(/r(t — r)) dr : Re/r < 0, 

< J —OO 


( 4 ) 

that guarantee the appropriate boundary conditions for 
t —> ±oo. The TS solution for any barrier motion is given 
by Eq. dSl). 


For anharmonic barriers (e ^ 0), the TS trajectory will 
be an unstable PO close to the barrier top, as in the har¬ 
monic case. Its period will typically coincide with the pe¬ 
riod T of the external driving. The anharmonic equations 
of motion are not amenable to an exact analytical 
solution, although approximate analytical methods have 
previously been employed42ii^ Instead we obtain the TS 
trajectory Fl = {x^(t), {t)) in phase space numerically 

as the periodic solution to the system of equations ©. A 
DS that is attached to Fl will be recrossing free. Phase 
space portraits of F^ are shown in Fig. [2] 


The barrier crossing rates for Eq. © were calculated 
by simulating ensembles of trajectories driven by an ex¬ 
ternal field of the form E{t) = asin(nt -|- (p). For sin¬ 
gle mode sinusoidal driving, the TS trajectory is a PO 
with period Physical units were set by normal¬ 

izing a and Wb to unity, making all other parameters 
dimensionless. Each trajectory was given an initial po¬ 
sition xq = —0.1 to the left of the instantaneous bar¬ 
rier top and Vo was sampled from two separate distribu¬ 
tions: (1) a Boltzmann distribution qb with ksT = 1, and 
(2) a uniform distribution gu (bounded over the region 
[pt — 1/2, Pt -1-1/2]). For each parameter set {0,7, e}, 
10®-10® trajectories were simulated. The normalized re¬ 
actant population PrII) is obtained from a histogram of 
those trajectories that are on the reactant side of the 
TS trajectory at time t. Assuming first order kinet¬ 
ics, the scaled logarithm of the normalized population, 
— In [PR{t) — Pr(oo)], should be linear in time. As illus¬ 
trated in Fig. [21 after transient trajectories have crossed, 
the decay of the logarithmic population is linear up to 
periodic modulation, and the first order assumption is 
confirmed. Periodic fluctuations are noticeable for small 
driving frequencies (fl ^ 2) and large anharmonicities 
due to effect of higher order terms in the asymptotic 
decay of Pnit). The slope of a least squares fit to the 
non-transient section of the data gives the reaction rates 
calculated from simulation k{. 


We now focus on analysis of the TS trajectory and 
the determination of reaction rates from its intrinsic sta¬ 
bility. With the bounded TS trajectory now defined, a 
relative coordinate system can be introduced. In relative 
coordinates 



t 


FIG. 3. Time dependence of the scaled logarithm of the reac¬ 
tant popnlation for = 2 and fl = 7 with vo sampled from gs ■ 
The slope of each dashed line is the barrier crossing rate fcf, 
corresponding to a respective e value. Parameters are 7 = 1, 
and (p = 0 


the equations of motion read 
Ax = Av, 

Av = —jAv — P'(Ax + x^(t)) + P'(x^(t)). ^ ^ 

The last term represents a time-dependent driving for 
the relative dynamics that does not depend on the cur¬ 
rent trajectory. It ensures that the relative equations of 
motion have a fixed point AF* at Ax = Av = 0, i.e., on 
the TS trajectory. 

The long-time decay rate of Pnit) is determined by 
the behavior of trajectories close to the stable manifold. 
Once a trajectory is sufficiently close to the TS trajectory, 
it can be described by a linearization of the equations of 
motion (O, 

Ax = Av, ^ ^ 

(7) 

Av = —jAv — a(t) Ax, 


where a(t) = [/"(x^t)). In the phase space vector coor¬ 
dinate AF = (Ax, Av) this linearization is given by 


AF = J(t) AF 

where 


J(t) = 


0 


( 8 ) 


Ax = x — x^(t), Av = v — v^(t), 


( 5 ) 


-I- 3e{xi{t) — E(t))‘^ 


-7 


(9) 
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is the Jacobian of Eq. (O about AF*. The linearity of 
Eq. ([5|) allows its solution to be expressed as 

Ar(t) = cr(t) Ar(T) (10) 

where the fundamental matrix solution cr{t) is a 2 x 2 
matrix that satisfies 

a- = J{t)cr, cr{0) = I, (11) 

where I is the identity matrix. 

The fundamental matrix for one period of AF^ is the 
monodromy matrix M = cr{T) whose eigenvalues TOu.s 
are called Eloquet multipliers. The Floquet exponents 
Mu.s = 1/Tln |mu s| give the rates by which nearby tra¬ 
jectories approach or recede from AF^i^ For a har¬ 
monic barrier, the multipliers are bonded according to 
0 < TOs < 1 < rUu giving rise to a positive Floquet expo¬ 
nent and a negative exponent /Tg. We will assume that 
this qualitative condition is also satisfied for the anhar- 
monic barriers; we neglect the possibility that for strong 
anharmonicities bifurcations of the TS trajectory might 
occur. 

Let t’u.s(O) be the eigenvectors of M. By Floquet’s 
theorem and the positivity of the Floquet multipliers, 
the vectors 


Vu,s{t)=e cr(t)t>u_g(0) (12) 


are periodic in time with period T. In the coordinate 
system defined by these vectors. 


AF(t) = Zu{t) Vu{t) + Zg(t) Vs{t), 

(13) 

the linearized equations of motion (|S]) read 


- /^U,s ^U,S5 

(14) 

with the solution 


Zn,sP) = 

(15) 


Therefore, the vectors Vu,s(i) determine the instanta¬ 
neous directions of the stable and unstable manifolds in 
the linear approximation. The actual stable and unstable 
manifolds are tangent to these directions at the TS tra¬ 
jectory. 

According to Eq. (ITSl) . the dynamics of Eq. (O is there¬ 
fore given by 

Ax{t) = Cu auit) -|- Cg ag(t) e'""*, (16) 


where Ou.s are the first components of the vectors Vu_g. 
They are periodic with period T. A trajectory with given 
initial conditions Cu and Cg will cross the moving divid¬ 
ing surface Ax = 0 at time t determined by 


Cu Ou(t) 


(17) 


If the initial condition Cg is fixed and a trajectory 
with a certain value of Cu crosses the moving DS at 



=1 
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FIG. 4. The barrier crossing rates for Hamiltonian (top) and 
dissipative (bottom) systems following Eq. ([51 . The numeri¬ 
cally calculated rates fcf for the distributions qs (circles) and 
qu (squares) are shown with units given on the left axes. The 
solid curves are the rates predicted by the difference in the 
Floquet exponents pu — Ps of the TS trajectory with units 
given on the right axes. 


time t, Eq. (II3 shows that a trajectory with initial value 
will cross at time t + T. Iteration then 
leads to the existence of trajectories with initial values 
that cross at time t + nT. 

Now consider an arbitrary ensemble of initial condi¬ 
tions with fixed x(0) on the reactant side and with a 
fixed value Cg < 0 small enough to be in the region of 
phase space where the linear approximation (|S]) is accu¬ 
rate. In this region the phase space density is constant up 
to linear corrections in the distance from the stable man¬ 
ifold and the number of trajectories that cross the DS in 
a given time interval is proportional to the width of the 
strip that contains these trajectories. From one period 
to the next this width decreases by a factor 
Thus, up to periodic modulation, the flux must decay 
by this same factor. The flux through the moving DS 
is the time derivative of the population, Fuit) = P{t), 
and thus the decay of PYi{t) is proportional to 
From this decay rate it follows that, fcf = ^u ~ Msj which 
states that the rate of barrier crossing is the difference 
in the Floquet exponents. Note that we have made no 
assumption for the energy distribution and thus this rate 
is independent of the ensemble of initial conditions. 

A comparison between the rates calculated from nu¬ 
merical simulation fcf, for both the Boltzmann qe and 
uniform qu distributions, and rates predicted by the Flo- 
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quet exponents /^u ~ Ms is shown in Fig. |4l For all values 
of the forcing frequency fl, dissipative parameter 7, and 
anharmonic strength e, the numerical rate is in agree¬ 
ment with rate predicted by stability analysis. This re¬ 
sult opens the possibility that when chemical reactions 
are forced by periodic external fields the reaction rates 
can be extracted from knowledge of the stability of the 
TS trajectory. The extension of TS trajectory stability 
analysis to aperiodically forced or thermally activated re¬ 
actions is a focus of our future research. 
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